Librairies et jeu de données

library(rmarkdown)
library(tidyverse) # pour  %>%
library(knitr) # kable
library(e1071) # skewness et kurtosis
library(ggplot2)
library(plotly) # graphes dynamiques
library(MASS)
library(fitdistrplus) #graphe Cullen & Frey

Ex: EPO .

Télécharger le jeu de données EPO modifié et nommer le mydata.

mydata=read.csv2("mydata2.csv",sep=',',row.names = 1)
paged_table(mydata)

Transformer la nature des variables nécessaires.

mydata$sexe<-as.factor(mydata$sexe)  
mydata$taille<-as.numeric(mydata$taille)   
mydata$Score<-as.factor(mydata$Score) 
mydata$niveau<-ordered(mydata$niveau,levels=c("L1_2","L3","M1_2"))
 str(mydata)
## 'data.frame':    38 obs. of  6 variables:
##  $ sexe             : Factor w/ 2 levels "F","H": 1 1 1 1 1 1 1 1 1 1 ...
##  $ age              : int  22 20 25 34 20 19 18 21 21 21 ...
##  $ taille           : num  1.7 1.66 1.65 1.62 1.66 1.64 1.61 1.63 1.63 1.65 ...
##  $ niveau           : Ord.factor w/ 3 levels "L1_2"<"L3"<"M1_2": 2 2 3 3 2 1 2 2 2 2 ...
##  $ scoreConnaissance: int  3 3 1 3 2 1 1 2 3 3 ...
##  $ Score            : Factor w/ 3 levels "=3","inf à 3",..: 1 1 2 1 2 2 2 2 1 1 ...

Pour alléger les écritures, procédons à des copies de nos variables

sexe<-mydata$sexe
age<-mydata$age
taille<-mydata$taille
niveau<-mydata$niveau
Score<-mydata$Score
scoreNum<-mydata$scoreConnaissance

1. Vocabulaires et premières définitions

  • Série statistique univariée : la suite \(x_1,\dots,x_n\).

  • Modalités: Valeurs distinctes prises par la série: \(m_1,\cdots,m_J\), \(J\leq n\).


Effectifs de la modalité \(m_j\): \[n_j\leq n \ \text{et} \ \sum_{j=1}^Jn_j=n.\]


Pour afficher le tableau des effectifs, on utilise la fonction table.

T_eff_niveau<-table(niveau) 
T_eff_niveau
## niveau
## L1_2   L3 M1_2 
##    6   28    4

Effectifs cumulés de la modalité \(m_j\) \[N_j=\sum_{k=1}^j n_k \ \text{et} \ N_J=n.\]


Pour afficher le tableau des effectifs cummulés, on utilise applique la fonctioncumsum au tableau des effectifs. P

cumsum(T_eff_niveau)
## L1_2   L3 M1_2 
##    6   34   38

Task 1

Task 1.A

Afficher la table des effectifs de la variable sexe et le nommer T_eff_sexe

Solution

T_eff_sexe<-table(sexe)
T_eff_sexe
## sexe
##  F  H 
## 23 15

Task 1.B

Afficher la table des effectifs cummulés de la variable sexe

Solution

cumsum(T_eff_sexe)
##  F  H 
## 23 38

Fréquence de la modalité \(m_j\): \[f_j=\frac{n_j}{n} \ \text{et} \ \sum_{j=1}^Jf_j = 1\]


prop.table(T_eff_niveau)
## niveau
##      L1_2        L3      M1_2 
## 0.1578947 0.7368421 0.1052632

Fréquence cumulée de la modalité \(m_j\): \[F_j=\sum_{k=1}^j f_k \ \text{et} \ F_J=1.\]


cumsum(prop.table(T_eff_niveau))
##      L1_2        L3      M1_2 
## 0.1578947 0.8947368 1.0000000

Task 2

Task 2.A

Afficher les Fréquences de la variable sexe

Solution

prop.table(T_eff_sexe)
## sexe
##         F         H 
## 0.6052632 0.3947368

Task 2.B

Afficher les Fréquences cummulées de la variable sexe .

Solution

cumsum(prop.table(T_eff_sexe))
##         F         H 
## 0.6052632 1.0000000

2. Cas des variables quantitatives continues : Choix des classes.

Lorsque le nombre de modalités trop grand (\(\Rightarrow\) pas de sens), il faut la traiter comme une variable continue et définir des classes.

table(taille) 
## taille
##  1.6 1.61 1.62 1.63 1.64 1.65 1.66 1.67 1.68 1.69  1.7 1.72 1.73 1.74 1.85 
##    1    1    2    3    4    3    4    2    3    2    5    4    2    1    1

A l’aide de la fonction cut( , breaks= ), on peut choisir le nombre de classes de même taille pour faire le tableau d’effectif.

T_eff_taille_5<-table(cut(taille,breaks=5))
T_eff_taille_5
## 
## (1.6,1.65] (1.65,1.7] (1.7,1.75] (1.75,1.8] (1.8,1.85] 
##         14         16          7          0          1

Task 3

Task 3.A

Afficher le tableau d’effectifs de la variable taille en considérant 4 classes et le nommer T_eff_taille_4.

Solution

T_eff_taille_4<-table(cut(taille,breaks=4))
T_eff_taille_4
## 
##  (1.6,1.66] (1.66,1.73] (1.73,1.79] (1.79,1.85] 
##          18          16           3           1

Task 3.B

Afficher les fréquences de la variable taille en considérant 4 classes.

Solution

prop.table(T_eff_taille_4)
## 
##  (1.6,1.66] (1.66,1.73] (1.73,1.79] (1.79,1.85] 
##  0.47368421  0.42105263  0.07894737  0.02631579

3. Principaux indicateurs

Les indicateurs décrivent synthétiquement la distribution des données.

⚠️ Attention !!

Un indicateur ne peut être pris seul comme description de la distribution mais doit être comparé aux autres. En outre, les paramètres de position et de dispersion doivent être mis en regard pour avoir une idée plus complète de la distribution.


3.1. Indicateurs de tendance centrale et de position

Moyenne empirique. Pour les variables quantitatives uniquement! \[\overline{x}= \frac{x_1+ \dots +x_n}{n}=\frac{1}{n} \sum_{i=1}^n x_i=\sum_{j=1}^Jf_jm_j.\]


mean(age)
## [1] 21.44737

Mode(s)/ Classe(s) modale(s) empirique(s). Modalité(s) la/les plus fréquente(s). Pas nécessairement unique!


Pour les variables catégorielles ou discrètes, le(s) mode(s) est la (les) modalité(s) la (les) plus fréquente(s):

c(which.max(table(niveau)),max(table(niveau)))
## L3    
##  2 28

Ici, le mode est unique et est la modalité numéro 2, soit la modalité L3 avec un effectif de 28

💡 Bon à savoir !

Pour les variables continues, on parle de classe(s) modale(s) Pas nécessairement unique!. Pour l’identifier, il faut caluler les fréquences renormalisées. Parmi les Classes \([b_{i-1},b_i[\), la (ou les) classe(s) modale(s) est (sont) celle(s) qui a (ont) la (les) fréquence(s) renormalisée(s) \(h_i\) la plus grande. \[h_i=f_i/(b_i-b_{i-1}).\] Si toutes les classes ont la même taille \((b_i-b_{i-1})\) alors il suffit de regarder les fréquences.


T_eff_taille_5<-table(cut(taille,breaks=5))
c(which.max(T_eff_taille_5),max(T_eff_taille_5))
## (1.65,1.7]            
##          2         16

Ici, la classe modale est unique et est la seconde classe, soit la classe (1.65,1.7] avec un effectif de 16

⚠️ Attention !

Si on change les nombre de classes, cela change la classe modale!


T_eff_taille_8<-table(cut(taille,breaks=8))
c(which.max(T_eff_taille_8),max(T_eff_taille_8))
## (1.63,1.66]             
##           2          11

Ici, la classe modale est unique et est la seconde classe, soit la classe (1.62,1.65] avec un effectif de 10.

Distribution unimodale : un seul maximum marqué dans le diagramme en bâtons ou l’histogramme.

Distribution multimodale : plusieurs maxima relatifs. Signe d’hétérogénéité de la population par rapport à la variable étudiée. (\(\rightarrow\) Existence de sous-populations).


Les quantiles empiriques Soit \(\alpha\in (0,1)\) une proportion que l’on se fixe. Le quantile d’ordre \(\alpha\), noté \(q_\alpha\) est tel que :(\(\alpha\%\) des \(x_i\leq q_{\alpha}\)) et (\(1-\alpha\%\) des \(x_i\geq q_{\alpha}\)).}


Variables quantitatives : quelques quantiles.

  • Si \(\alpha=0.25\), alors \(q_1=q_{0.25}\) est appellé 1er quartile (quantile(.,prob=0.25))

  • Si \(\alpha=0.75\), alors \(q_3=q_{0.75}\)est appellé 3ème quartile (quantile(.,prob=0.75))

  • Si \(\alpha=0.5\), alors \(me=q_2=q_{0.5}\)est appellé 2d quartile (median(.))

La fonction quantile()permet d’afficher les 3!


quantile(age)
##   0%  25%  50%  75% 100% 
##   18   20   21   22   34

📌 Pour les variables qualitatives ordonnées, il est possible de calculer des quantiles. Le plus simple est de convertir temporairement la variable en variable numeric, puis de calculer le quartile à l’aide de la commande level()


levels(niveau)[median(as.numeric(niveau))]
## [1] "L3"

💡 Bon à savoir !

  • Médiane \(vs\) Moyenne

    • Médiane souvent plus intéressante que moyenne. Individu médian existe et non nécessairement l’individu moyen.
    • Médiane est peu sensible aux valeurs extrêmes.


  • Position des indicateurs pour une distribution unimodale

    • Distribution symétrique: \(m_e=\overline x\) et la classe modale contient la médiane et la moyenne.
    • Distribution asymétrique: \(m_e\) souvent située entre la classe modale et \(\overline x\).


Task 4

Task 4.A.

Calculer la taille moyenne et la médiane.

Solution

paste0('mean = ',mean(taille))
## [1] "mean = 1.67578947368421"
paste0('median = ', median(taille))
## [1] "median = 1.67"

Task 4.B.

Donner la classe modale pour un découpage en 4 classes de la variables taille.

Solution

T_eff_taille_4<-table(cut(taille,breaks=4))
c(which.max(T_eff_taille_4),max(T_eff_taille_4))
## (1.6,1.66]            
##          1         18

Ici, la classe modale est unique et est la première classe, soit la classe (1.6,1.66] avec un effectif de 18.

Task 4.C.

Donner Calculer le mode de la variable Score.

Solution

table(Score)
## Score
##      =3 inf à 3 sup à 3 
##      15      14       9
c(which.max(table(Score)),max(table(Score)))
## =3    
##  1 15

Le mode est la modalité =3avec un effectif de 15.

Task 4.D.

Afficher le min, le max et les 3 quartiles pour la variable scoreConnaissance.

Solution

quantile(mydata$scoreConnaissance)
##   0%  25%  50%  75% 100% 
##    0    2    3    3    6

Task 4.E.

Pour la variable qualitative ordonnée niveau, calculer le 1er quartile.

Solution

levels(niveau)[quantile(as.numeric(niveau),prob=0.25)]
## [1] "L3"

La commande summary() affiche un résumé des données qui diffère selon leur nature.

summary(mydata)
##  sexe        age            taille       niveau   scoreConnaissance     Score   
##  F:23   Min.   :18.00   Min.   :1.600   L1_2: 6   Min.   :0.000     =3     :15  
##  H:15   1st Qu.:20.00   1st Qu.:1.640   L3  :28   1st Qu.:2.000     inf à 3:14  
##         Median :21.00   Median :1.670   M1_2: 4   Median :3.000     sup à 3: 9  
##         Mean   :21.45   Mean   :1.676             Mean   :2.763                 
##         3rd Qu.:22.00   3rd Qu.:1.700             3rd Qu.:3.000                 
##         Max.   :34.00   Max.   :1.850             Max.   :6.000

3.2. Indicateurs de dispersion

Les indicateurs de dispersion sont uniquement définit pour des variables quantitatives. Ils informent sur la variabilité des observations autour d’une valeur centrale. Ils sont toujours positif et d’autant plus grands que la variable est dispersée (grande variabilité autour d’une caractéristique de tendance centrale).

Variance empirique, écart typeempirique.

  • Variance empirique est la moyenne des carrés des différences entre \(x_i\) et la moyenne arithmétique \(\overline{x}\) \[s_X^2 = \frac{1}{n}\sum_{i=1}^n (x_i-\overline{x})^2=\overline{x^2}-(\overline{x})^2.\]

  • Écart type empirique: \[s_X = \sqrt{s_X^2}.\]


En général, les logiciels de statistiques (R inclus) fournissent la variance (et écart type) corrigée \[(s_X^c)^2=\frac{n}{n-1}s_X^2.\]

rbind(paste0(c("variance corrigée: ", var(age))),
paste0(c("variance : ", mean(age^2)-(mean(age))^2)),
paste0(c("Écart type corrigé: ", sd(age))),
paste0(c("Écart type : ", sqrt(mean(age^2)-(mean(age))^2))))
##      [,1]                   [,2]              
## [1,] "variance corrigée: "  "7.01066856330014"
## [2,] "variance : "          "6.82617728531864"
## [3,] "Écart type corrigé: " "2.64776671240125"
## [4,] "Écart type : "        "2.61269540615025"

Intervalle interquartile, étendue

  • Etendue: différence entre la plus petite et la plus grande valeur de la série (Dépend uniquement des valeurs extrêmes) \[E = x_{(n)}-x_{(1)}=\max-\min.\]

  • Distance inter-quartiles :

\[IQ = q_{0.75}-q_{0.25}=q_3-q_1.\]


rbind(paste0(c("Etendu : ", max(age)-min(age))),
paste0(c("IQ : ",quantile(age,prob=0.75)-quantile(age,prob=0.25) )))
##      [,1]        [,2]
## [1,] "Etendu : " "16"
## [2,] "IQ : "     "2"

3.3. Indicateurs de forme

Les indicateurs de forme fournissent une information sur la forme de la distribution par des mesures de son asymétrie et de son applatissement.

  • Asymétrie par rapport à la moyenne : Coefficient d’asymétrie de Fisher.

  • Aplatissement d’une distribution (importance des “queues” de la distribution) : Coefficient d’aplatissement de Fisher.

Coefficient d’asymétrie de Fisher (Skewness)

\[\gamma_1=\frac 1n\sum_{i=1}^n\left(\frac{x_i-\overline{x}}{s_X}\right)^3\]

  • \(\gamma_1=0\) : distribution (unimodale) symétrique par rapport à sa moyenne.

  • \(\gamma_1<0\) : distribution unimodale étalée vers la gauche (queue de distribution à gauche).

  • \(\gamma_1>0\) : distribution unimodale étalée vers la droite (queue de distribution à droite).

📌 \(\gamma_1\) compare les moments de la distribution observée avec ceux de la distribution gaussienne (en “cloche”, symétrique) pour laquelle \(\gamma_1=0\).


library(e1071)
skewness(age)
## [1] 2.698495

\(\gamma_1>0\) donc distribution unimodale approximativement étalée vers la droite (queue de distribution à droite).

Coefficient d’aplatissement de Fisher (Kurtosis)

\[\gamma_2=\frac 1n\sum_{i=1}^n\left(\frac{x_i-\overline{x}}{s_X}\right)^4-3\]

  • \(\gamma_2>0\) : la distribution observée est moins aplatie qu’une distribution gaussienne (de même moyenne et de même variance).

  • \(\gamma_2<0\) : la distribution observée est plus aplatie qu’une distribution gaussienne (de même moyenne et de même variance).

  • \(\gamma_2=0\) : distribution gaussienne centrée réduite.


kurtosis(age)
## [1] 10.69388

\(\gamma_2>0\) donc distribution moins aplatie qu’une distribution gaussienne (de même moyenne et de même variance).

4. Représentations graphiques

  • Diagramme en barres : variable qualitative.

  • Diagramme en secteurs : variable qualitative.

  • Diagramme en batons : variable quantitative discrète.

  • Histogramme : variable quantitative continue.

  • Boxplot : variable quantitative discrète et continue.

📌 La commande plot() laisse le choix de la représentation à R.

Diagramme en barres
La hauteur de chaque barre est proportionnelle à l’effectif (ou à la fréquence). La largeur de chaque barre est arbitrairement choisie.

library(plotly)
plot_ly(x = names(table(mydata$niveau)), y = table(mydata$niveau), type = 'bar', marker = list(color = "purple")) %>%
  layout(title = 'Niveau: Diagramme en barres',
         xaxis = list(title = 'Niveau'),
         yaxis = list(title = 'Effectifs'))

Diagramme en secteurs
Chaque secteur angulaire est proportionnel à l’effectif (ou à la fréquence).

# Choisir la couleur des barres

couleur <- c("#FF5733", "#33FFB1", "#3379FF") 
couleur1 <- c("#890090", "#402023", "#59F9FF") 
plot_ly(labels =  names(table(mydata$niveau)), values = table(mydata$niveau), type = 'pie', marker = list(colors = couleur1)) %>%
  layout(title = 'Niveau: Diagramme en secteurs')

Diagramme en bâtons

Les batônnets (traits) ont pour abscisse la modalité et sont de hauteur proportionnelle à l’effectif (ou à la fréquence).

plot_ly(x = names(table(mydata$age)), y = table(mydata$age)) %>%
  layout(title = 'Age: Diagramme en bâtons',
         xaxis = list(title = 'Age'),
         yaxis = list(title = 'Effectifs'))%>%
  add_bars(x= names(table(mydata$age)),y= table(mydata$age),width = 0.01)

Histogramme
Dans un histogramme, les rectangles ont pour bornes les bornes de chaque classe et la hauteur est soit l’effectif soit la fréquence.

plot1<-plot_ly(x = mydata$taille, name = "Effectifs",
        type = 'histogram', 
        nbinsx = 10,  
        marker = list(color = "purple")) %>%
  layout(title = 'Niveau :Histogramme',
         xaxis = list(title = 'Taille'),
         yaxis = list(title = 'Effectif')) %>% 
  layout(
  barmode="overlay",
  bargap=0.1)



plot2<-plot_ly(mydata)%>% add_histogram(x = taille, name = "Frequences",
                             histnorm = "probability", 
                             marker = list(color = "cyan")) %>%
  layout(title = 'Niveau: Histogramme',
         xaxis = list(title = 'Taille'),
         yaxis = list(title = 'Fréquences')) %>% 
  layout(
  barmode="overlay",
  bargap=0.1)


# Affichage des histogrammes côte à côte
subplot(plot1, plot2)

Boxplot/Boîte de distribution/boîte à moustache
Résume quelques caractéristiques de position et de dispersion du caractére étudié (médiane, quartiles, minimum, maximum ou déciles).


#plot_ly(mydata)%>% add_boxplot(y = taille, marker = list(color = "cyan")) %>%
#  layout(title = 'Taille: Boxplot',
#         xaxis = list(title = 'Taille'))


plot_ly(y = mydata$taille, type = "box",name = "", marker = list(color = "pink"))%>%
 layout(title = 'Taille: Boxplot',
     yaxis = list(title = 'Taille'))

📌 A noter!

  • Si La boxplot est raccourcie au min \(x_{(1)}\) et au max \(x_{(n)}\) lorsqu’aucune observation n’arrive au delà des moustaches.

  • Faire un examen attentif des individus “hors norme” : valeur abérrante (erreur saisie…) ou extrême (éloignée du centre) , pas d’attitude systématique.

Task 5

Task 5.A.

Proposer une representation graphiquede la variable sexe.

Solution

plot_ly(x = names(table(mydata$sexe)), y = table(mydata$sexe), type = 'bar', marker = list(color = "purple")) %>%
  layout(title = 'sexe: Diagramme en barres',
         xaxis = list(title = 'Niveau'),
         yaxis = list(title = 'Effectifs'))

Task 5.B.

Proposer une representation graphique de la variable Score.

Solution

couleur <- c("#FF5733", "#33FFB1", "#3379FF") 
plot_ly(labels =  names(table(mydata$Score)), values = table(mydata$Score), type = 'pie', marker = list(colors = couleur)) %>%
  layout(title = 'Score: Diagramme en secteurs')

Task 5.C.

Proposer une representation graphique de la variable scoreConnaissance.

Solution

plot_ly(x = names(table(mydata$scoreConnaissance)), y = table(mydata$scoreConnaissance)) %>%
  layout(title = ' Diagramme en bâtons',
         xaxis = list(title = 'scoreConnaissance'),
         yaxis = list(title = 'Effectifs'))%>%
  add_bars(x= names(table(mydata$age)),y= table(mydata$age),width = 0.01)

Task 5.D.

Proposer une representation graphique de la variable Taille.

Solution

plot_ly(x = mydata$taille, name = "Effectifs",
        type = 'histogram', 
        nbinsx = 7,  
        marker = list(color = "purple")) %>%
  layout(title = 'Niveau :Histogramme',
         xaxis = list(title = 'Taille'),
         yaxis = list(title = 'Effectif')) %>% 
  layout(
  barmode="overlay",
  bargap=0.1)

Task 5.E.

Proposer une representation graphique de la variable age.

Solution

plot_ly(y = mydata$age, type = "box",name = "", marker = list(color = "pink"))%>%
 layout(title = 'AGE: Boxplot',
     yaxis = list(title = 'Taille'))

5. Loi empirique/ adéquation à une loi

Fonction de répartition empirique
Les fréquences cumulées sont représentées par une fonction en escalier.


library(ggplot2)
library(plotly)


# Création du graphique ECDF avec ggplot2
gg_ecdf <-
  ggplot(mydata, aes(x = age)) +
  stat_ecdf(geom = "step",color='blue') +
  labs(title = "Fonction de répartition empirique",
       x = "Âge",
       y = "F_n") +theme_light() 
# Affichage du graphique interactif 
ggplotly(gg_ecdf)

Diagramme Quantile-Quantile (QQ-plot)
Outil graphique permettant de visualiser rapidement l’adéquation de la distribution d’une série numérique à la distribution de référence \(F_0\). Dans ce graphe, on reporte sur l’axe des ordonnées les fractiles correspondant à la distribution observée et sur l’axe des abscisses ceux correspondant à la distribution théorique.




library(ggplot2)

library(ggplot2)

n <- 1000
m <- 5
sigma <- 2
x_lim <- c(-m, m)
y_lim <- x_lim

# Exemple de données
Normal_data <- data.frame(valeurs = rnorm(n, 0, 1))
sigma_data <- data.frame(valeurs = rnorm(n, 0, sigma))
m_data <- data.frame(valeurs = rnorm(n, m, 1))
m_sigma_data <- data.frame(valeurs = rnorm(n, m, sigma))
uni_data<-data.frame(valeurs =-runif(n, min = 0, max = 5))

# Création des graphiques QQ plot avec ggplot2 pour chaque jeu de données
ggplot() +
  stat_qq(data = Normal_data, aes(sample = valeurs, color = "Normal"), distribution = qnorm, dparams = list(mean = 0, sd = 1), size = 0.1) +
  stat_qq(data = sigma_data, aes(sample = valeurs, color = "Sigma=2"), distribution = qnorm, dparams = list(mean = 0, sd = 1), size = 0.1) +
  stat_qq(data = m_data, aes(sample = valeurs, color = "Mean=5"), distribution = qnorm, dparams = list(mean = 0, sd = 1), size = 0.1) +
  stat_qq(data = m_sigma_data, aes(sample = valeurs, color = "Mean=5,Sigma=2"), distribution = qnorm, dparams = list(mean = 0, sd = 1), size = 0.1) +
  stat_qq(data = uni_data, aes(sample = valeurs, color = "Uniform[0,5]"), distribution = qnorm, dparams = list(mean = 0, sd = 1), size = 0.1) +
  geom_abline(intercept = 0, slope = 1, color = "red", linetype = "dashed") + 
  labs(title = "QQ plot - Distribution normale",
       x = "Quantiles théoriques",
       y = "Quantiles observés",
       color = "Données") +  # Légende pour les couleurs
  coord_cartesian(xlim = x_lim, ylim = y_lim) +
  scale_color_manual(values = c("Normal" = "blue", "Sigma=2" = "green", "Mean=5" ="orange","Mean=5,Sigma=2"="pink", "Uniform[0,5]"="cyan")) 

💡 Bon à savoir !

  • Points alignés sur 1ère bissectrice \(\Rightarrow\) \(F_0\) est adaptée (courbe blue).

  • Points alignés sur une droite parallèle à la 1ère bissectrice \(\Rightarrow\) erreur sur les paramètres de position de \(F_0\) (courbe orange).

  • Points alignés sur une droite passant par l’origine mais inclinée par rapport à la 1ère bissectrice \(\Rightarrow\) erreur sur les paramètres de dispersion de \(F_0\) (courbe green).

  • Points alignés sur une droite ne passant pas par l’origine et inclinée par rapport à la 1ère bissectrice \(\Rightarrow\) erreur sur les paramètres de dispersion et de position de \(F_0\) (courbe pink).

  • Points non alignés sur une droite \(\Rightarrow\) \(F_0\) non adaptée (courbe cyan).


Prenons pour la suite le jeu de données carsde la librarie MASS qui contient

50 Données et 2 variables : speed :la vitesse en miles par heure et dist la distance nécessaire à l’arrêt d’un véhicule.

library(rmarkdown)
library(MASS)
data_cars <- as.data.frame(cars)
paged_table(data_cars)

Histogramme avec la fonction hist() de base pour aller plus vite

hist(data_cars$dist, freq = FALSE, col = "purple", main = "", xlab = "Distance", ylab = "Densité")

A quelle loi pensez-vous ? Une gaussienne?

# Fonction de densité de la loi normale
densite_normale <- function(x) dnorm(x, mean(data_cars$dist), sd(data_cars$dist))

# Générer des données pour la courbe de densité
x <- seq(min(data_cars$dist), max(data_cars$dist), length.out = 1000)
y <- densite_normale(x)

# Tracer la courbe de densité avec lines()
hist(data_cars$dist, freq = FALSE, col = "purple", main = "", xlab = "Distance", ylab = "Densité")
lines(x, y, col = "red", lwd = 2)

Pas génial !

Graphe de Cullen et Frey

Permet de trouver la meilleure loi usuelle parmi différentes lois usuelles à l’aide de la fonction descdist de la librairie fitdistrplus.

Interprétation graphe de cullen & Frey

  • Le point bleu du graphe représente notre échantillon.

  • On regarde sur le graphe : De quelle loi, ce point est-il le plus proche ?


library(fitdistrplus)
distance=data_cars$dist
descdist(distance, boot=100)

## summary statistics
## ------
## min:  2   max:  120 
## median:  36 
## mean:  42.98 
## estimated sd:  25.76938 
## estimated skewness:  0.806895 
## estimated kurtosis:  3.405053

Ici la loi gamma convient/explique (“fit”) le mieux les données! Vérifions

plot(fitdist(distance, "gamma"))

🙂 Evaluation graphique de l’heuristique : Loi Gamma “fit” bien ! A titre de comparaison, évaluons graphiquement l’adéquation avec la seconde loi la “plus” convenable: la loi normale.

plot(fitdist(distance, "norm"))

La Loi normale “fit” moins bien !

Estimation des paramètres de la loi gamma \(\mathcal G(k,\theta)\)

Estimateur du maximum de vraisemblance : \((\widehat k,\widehat\theta)\) estimation de \((k,\theta)=(`shape,rate')\) est calculé par R. Pour afficher des intervalles de confiance à 95%, on utilise la commande confit()


fitdist(distance, "gamma")
## Fitting of the distribution ' gamma ' by maximum likelihood 
## Parameters:
##         estimate Std. Error
## shape 2.37677270 0.44562696
## rate  0.05529455 0.01153582
confint(fitdist(distance, "gamma")) 
##            2.5 %     97.5 %
## shape 1.50335990 3.25018550
## rate  0.03268476 0.07790435

Task 6

Task 6.A.

Quelle loi convient à la variable speed?

Solution

speed=data_cars$speed
descdist(distance, boot=100)

## summary statistics
## ------
## min:  2   max:  120 
## median:  36 
## mean:  42.98 
## estimated sd:  25.76938 
## estimated skewness:  0.806895 
## estimated kurtosis:  3.405053

La loi normale semble bien “fitter”

plot(fitdist(speed, "norm"))

Task 6.B.

Les données speed sont-elles issues d’une loi normale ? Proposer un test

Solution

Test de Shapiro \(H_0\): Loi normale \(v.s.\) \(H_1\): Pas une loi normale.

shapiro.test(speed)
## 
##  Shapiro-Wilk normality test
## 
## data:  speed
## W = 0.97765, p-value = 0.4576

\(\Rightarrow\) Petite \(p\)-value\(>5\%\), on ne rejette pas \(H_0\) donc oui, les données speed sont issues d’une loi normale!